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ABSTRACT 

The solar tachocline at the bottom of the convection zone is an important 
region for the dynamics of the Sun and the solar dynamo. In this region, the 
sound speed inferred by global helioseismology exhibits a bump of approximately 
0.4% relative to the standard solar model. Global helioseismology does not pro- 
vide any information on possible latitudinal variations or asymmetries between 
the Northern and Southern hemisphere. Here, we develop a time-distance helio- 
seismology technique, including surface- and deep-focusing measurement schemes 
and a combination of both, for two-dimensional tomographic imaging of the solar 
tachocline that infers radial and latitudinal variations in the sound speed. We 
test the technique using artificial solar oscillation data obtained from numerical 
simulations. The technique successfully recovers major features of the simplified 
tachocline models. The technique is then applied to SOHO/MDI medium-^ data 
and provides for the first time a full two-dimensional sound-speed perturbation 
image of the solar tachocline. The one-dimensional radial profile obtained by 
latitudinal averaging of the image is in good agreement with the previous global 
helioseismology result. It is found that the amplitude of the sound-speed pertur- 
bation at the tachocline varies with latitude, but it is not clear whether this is 
in part or fully an effect of instrumental distortion. Our initial results demon- 
strate that time-distance helioseismology can be used to probe the deep interior 
structure of the Sun, including the solar tachocline. 
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Introduction 



The solar tachocline, first defined by lSpiegel &: Zahnl (119921 ). is a region located between 
the solar convection zone and the radiative zone. It is a transition zone where the latitudinal 
differential rotation of the convection zone change s to a rigid rotation of the radiative zone, a s 



found by numerous helioseismology studies (e.g., iThompson et al.lll996t ISchou et al.lll998l ). 
In this region the inferred sound-speed per turbation profile exhibits the large st difference 
relative to a standard solar model - Model S (IChristensen-Dalsgaard et al.lll996r). as re vealed 



by a number of helioseismology studies (e.g.. lGough et al 



1996 



Kosovichev et al.lll997l ). The 



solar tachocline is also widely believed to be the location where the solar dynamo operates, 
since the sharp rotational gradient observed in this area can build up a strong toroida l 



field 



,e.g. 



Parker! Il993l ; ICharbonneau fc MacGregorl 119971 ; iDikpati Charbonneaul Il999l ) . 



Therefore, it is important to measure the structure, dynamics, and cyclic variations of the 
tachocline. Measurements in this zone will provide a better understanding of the transition 
from the solar convection zone to the radiative zone, the mechanisms of solar magnetic field 
generation and amplification, as well as the periodicity and evolution of the solar activity 
cycles. 

However, despite the importance of the tachocline, little is known about this area from 
observations besides the steep rotational shear and the one- dimensional sound-speed per- 
turbation bump. This region is difficult to probe by helioseismology because it is quite 
thin, ~ 0.05 — O.O9i? ( IKosovichevI Il996bl ; lBasulll997l). and is located deep in the Sun, ap- 
proximately O.29i? below the photosphere (IChristensen-Dalsgaard et al.lll99ll ). Nearly all 
observational results related to the tachocline area were obtained by global helioseismology 
analysis of normal mode frequencies, which is not capable of resolving hemispheric asym- 
metries of the differential rotation and sound-speed structures. In this paper, we explore 
the possibility of mapping a two-dimensional image of the ta chocline, both latit ude and 
depth dependent, by the use of time-distance helioseismology (IDuvall et al.lll993l ). a local 
helioseismology technique. 

Because helioseismology infers solar interior properties that are undetectable by other 
independent observations, i t is often difficult t o asse ss th e quality an d accu racy of the in- 
ferred properties. Recently, iGeorgobiani et al.l (120071 ) and IZhao et al.l (120071 ) evaluated the 
performance of the time-distance technique for inferring sma ll-scale subsurface fl ow fields by 
analyzing data from numerically si mulated solar c onvec tion (jBensen et al. 2006 ). and found 
satisfactory results. More recen tly, lHartlep et al.l (120081 ) validated the time-distance far-side 
imaging technique (jZhaol 120071 ) by using numerical simulations of the global acoustic wave- 
fields of the whole Sun. In this paper, we use a similar approach: develop time-distance 
measurement and inversion procedures for imaging sound-speed variations in the tachocline, 
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test these procedures on numerical simulation data, and finally apply the procedures to 
actual observations. 

To image the solar tachocline, we develop two time-dis tance measurement procedures 
based on surface- and deep-focusing schemes (jDuvalll 120031 ) . and corresponding inversion 
codes. We also invert combined measurements of these two schemes to improve the resolution. 
A description of the measurement and inversion techniques is given in Sec. 2. In Sec. 3 
we describe the application of our method to numerical simulation data, and assess the 
performance of our method. In Sec. 4, we present initial results of the tomographic imaging 
of the solar tachocline using SOHO/MDl data fjScherrer et al.lll995l ). A brief discussion 
follows in Sec. 5. 



2. Techniques 

2.1. Surface- and Deep-Focusing Measurements 

Time-distance helioseismology measures travel times of acoustic waves that propagate 
from one solar surface location to another along curved paths through the solar interior. 
These measurements are done by computing the cross-correlation between acoustic signals 
observed at these two separate locations. In practice, in order to improve the signal to noise 
ratio, the cross-correlations are calculated between the acoustic signal observed at a given 
point and the signal averaged over its surrounding annulus or part of the annulus. When all 
ray paths converge at one surface location in such measurements, as shown in Fig. [Tfe,, this 
measurement scheme is often called surface-focusing. When all ray paths pass through one 
point located in the solar interior, this scheme is often referred to as deep-focusing (Fig. [Tb). 

The measurement procedure is as follows: A datacube of a time sequence of surface 
velocity images, either numerically simulated or observed, is remapped into heliographic co- 
ordinates (using Postel's projection) with the center of the projected map located at the solar 
disk center and a pixel size of 0?6 (hereafter, degree means heliographic degree), consistent 
with the MDI medium-^ data sampling resolution. For the surface-focusing scheme every 
pixel in the dataset serves as a central point, and for each central point the acoustic signals 
are averaged over a set of annuli around this point. Each annulus corresponds to a particular 
range of travel distances of acoustic waves. The radius of the annuli ranges from 6?0 to 84?0 
with a step of 0?6, hence a total of 131 annuli (travel distances) are used. This number is 
significantly greater than the number of annuli used in most previous time-distance helio- 
seismology measurements. At each location, the cross-correlations are computed between 
the signals at the central point and the signals from each annulus. To further improve the 
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signal-to-noise ratio, we average the cross-correlations in overlapping latitudinal zones 5?0 
wide, from —60° to 60° with a step of 0?6. Then, we use the Gabor wavelet fitting procedure 



( IKosovichev &: Duvalllll997f ) to derive the mean acoustic phase travel times for all annuli at 



each latitude. 

For the deep-focusing scheme (Fig. [TJd), every pixel also serves as a central point, but 
the annuli are divided into four quadrants in the North, East, South, and West directions. 
The acoustic signals are averaged in each quadrant, and cross- correlations are calculated 
between opposite quadrants. Then, these two cross-correlations are averaged together. Note 
that for some locations, not both cross-correlations can be obtained, since some quadrant 
annuli may be located beyond the solar limbs or poles. In these measurements, the distance 
between the opposite quadrants varies from 6?0 to 84?0 with a step of 1?2; hence, a total 
of 66 travel distances is used. Similar to the surface-focusing scheme, the cross-correlations 
are averaged in longitudinal belts, and the Gabor wavelet fitting is performed to derive the 
acoustic travel time for each annulus and each latitudinal belt. 

It is important to note that unlike in most other time-distance helioseismology studies, 
the phase-speed filtering is not applied in this study, because for the large distances used 
i n this study , such filtering may introduce systematic shifts that are not fully understood 



( IDuvalll 120031 ) . However, as usual, solar convection and /-modes are removed by filtering. 



2.2. Inversion Algorithm 



After the acoustic travel times are measured an inversion procedure is applied to infer 
the sound-speed perturbations inside the Sun. The travel-time sensitivity kernels in sph erical 
coordinates are calculated using the ray-path approximation (jD'Silva fc Duvalllll995l ). The 
sensitivity kernels are computed separately for the two different measurement schemes. First, 
three-dimensional kernels are computed for all different annuli, and the kernels are then 
averaged along the longitudinal direction to make the kernels depend only on latitude and 
radius. 



The inversions are performed by solving linearized equati ons for the travel time varia- 
tions in the sense of least squares (IKosovichev fc Duvalllll997t ). Generally, two different in- 
version algorithms were used in the past in time-distan ce helioseismology, LSQR algorithm 



( IKosovichev! Il996al ; IZhao. Kosovichev. fc Duvall 



also known as MCD (e.g. 
employ the MCD technique 



Jacobsen et al. 



199 



9; 



200 if) and Multi-Cha nnel Deconvolution 



Couvidat et al.ll2005l ). In this paper, we 
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Results from Numerical Simulation Data 



3.1. Model and Numerical Simulation 



A numerical simulation code has been developed to sim ulate the three-dimensiona l 
acoustic wavefields of the global Sun in spherical coordinates ( lHartlep &: Mansourl 120081 ) . 
It solves the linearized Euler equations describing wave propagation excited by random 
sources near the surface. In the present study, the simulations include oscillations of spher- 
ical degree, £, in the range from to 170. This is sufficient for most global scale studies. 
Comparing with observations obtained from the MDI medium-^ program, this simulation 
code reproduces accurately the solar oscillation power spectrum {k — u diagram) except at 
high frequencies. It also reproduces the acous tic travel times measured by the time-distance 



technique (see Fig. 1 in lHartlep et al.l 120081 ). Global wavefields generated with this sim- 



ulation co de have been succes sfully used for validating the time-distance far-side imaging 
technique (lHartlep et al.l 120081 ) . 



The sound-speed perturbation near the solar tachocline ob served by global helioseismol 



ogy i s measured relative to a standard solar model, Model S (jChristensen-Dalsgaard et al. 



19961 ). Therefore, to infer the sound-speed structure at the tachocline from the measured 
acoustic travel times, it is necessary to subtract the corresponding reference travel times 
measured from a simulation for the unperturbed Model S. A latitudinally averaged refer- 
ence time-distance curve obtained for the surface-focusing scheme is shown in Fig. [2j The 
time-distance relationship for the deep-focusing scheme is similar. 

For testing, we use an artificial model of the sound-speed perturbations in the tachocline 
shown in Fig. [3j The perturbed structure has a two-dimensional Gaussian shape in each 
hemisphere, but is asymmetric relative to the equator. The maximum amplitudes of sound- 
speed perturbation, 0.7% in the Northern hemisphere and 0.5% in the Sou thern hemisphere. 



are n ot far apart from the values derived from global helioseismology (e.g.. lKosovichev et al. 
19971 ). The perturbation in each hemisphere is centered at 0.70i? Q with a FWHM of 0.082i? Q 
(the variance of the Gaussian is 0.035i?©) in the radial direction, and centered at latitude of 
30° with a FWHM of 35° (variance of 15°). The time series from the numerical simulations 
of the acoustic wavefields for both unperturbed and perturbed models are 1024 minutes long. 



3.2. Surface-Focusing Results 



For the surface-focusing scheme the measured acoustic travel times after subtraction of 
the reference travel times are shown in Fig. as a function of the annulus radius (travel 
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distance) and latitude. One can easily see that the travel times are apparently smaller 
on the right-hand side of the image when the radius is larger than 46°. This can be better 
identified in the ID curves obtained after latitudinal averaging in both hemispheres, as shown 
in Fig. Hb. 

The inversion results for the surface-focusing scheme (Fig. Hb,d) show that the model 
sound-speed perturbation in the tachocline area is recovered, although it does not perfectly 
match the two Gaussian shaped structures in both hemispheres. The location and amplitude 
of the perturbation match well. However, it is obvious from the sound-speed image (Fig. [4b ) 
and the latitudinally averaged ID sound-speed perturbation curves from both hemispheres 
(Fig. HJi), that the inverted sound-speed perturbation is wider than the prescribed model, 
especially into the deeper interior. A negative sidelobe close to the surface can also be seen. 

It is not immediately clear why the inverted sound-speed results are not well localized 
in the deeper area. Inversion averaging kernels, which characterize the spatial resolution, 
are displayed in Fig. [5] for three different targeted depths inside the Sun. It can be seen 
that when the target is located near the bottom of the convection zone the averaging kernel 
is more spread in latitude, and that the kernel is more spread in the radial direction when 
the target location is shallower. Nevertheless, the averaging kernels are well localized and 
should not be the reason for the spreading of the inverted sound-speed perturbation. We 
suspect that the poor localization with depth below the tachocline may be caused by the 
use of ray-path sensitivity kernels and by relatively high realization noise in the simulation 
data, which is significantly higher than in the real data since the simulation time series are 
shorter. The inversion m ay be improved by the use of Born-approximation kernels (e.g., 



Birch fc Kosovichevll2000l ). longer time series, and more realistic wave excitation models. 



3.3. Deep-Focusing Results 

The measurement and inversion results for the deep-focusing scheme are presented in 
Fig. El Similar to the surface-focusing measurements, a travel-time drop can be easily iden- 
tified at an annulus radius of ~ 45° from both the 2D travel-time image (Fig. [6^) and the 
ID latitudinally averaged curves (Fig. [6b). It is noteworthy that at high latitudes and large 
annulus radii, the acoustic travel-time measurements are often unreliable because too few 
pixels are available for averaging. The measurements in these areas are discarded to make 
sure that the inversion results are not contaminated by bad fitting points. 

The inversion results for the deep-focusing scheme are shown in Fig. [6b, d. Similar 
to the surface-focusing inversion results, the sound-speed bump is successfully recovered, 
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although the shape does not perfectly match the model. Again, the inferred sound-speed 
perturbation is spread into areas below the tachocline region. Despite that the averaging 
kernels for the deep-focusing inversions (Fig. [7]) are well localized, even better than for the 
surface-focusing inversions, the noise level of the inverted results is higher. Of course, there 
is a trade-off between localization and noise amplification in the inversion results. In the 
model calculations, the trade-off parameter is chosen to get good localization in the tachocline 
region. We did not try to optimize the inversion procedure for the deeper interior. The same 
regularization parameters are used for inversion of the real data. 



3.4. Results from a Combination of Surface- and Deep-Focusing 

The surface- and deep-focusing measurement schemes use the same datasets but are 
different analysis procedures, and have different sensitivities and signal-to-noise ratios at 
different depths. Thus, it is useful to combine the acoustic travel time measurements from 
both schemes, and make one inversion to improve the results. It is expected that such joint 
inversion would improve localization of the inversion results and enhance the signal-to-noise 
ratio. It is not difficult to revise the inversion code to solve two sets of linear integral equations 
instead of one using a proper combination of the surface- and deep-focusing sensitivity kernels 
and measured travel times. 

The inversion results are shown in Fig. [HJ As it can be seen, both the 2D and ID results 
look similar to an average of the separate inversion results of the surface- and deep-focusing 
schemes. Generally, the inversion from the combined measurements is closer to the model 
than each separate inversion. The averaging kernels from this inversion also look similar to 
the average of the averaging kernels of these measurement schemes, and are not shown here. 



4. Results from SOHO/MDI Observations 
4.1. Data Processing 



The medium-^ program of the MDI instrument onboard SOHO (jKosovichev et al.lll997l ) 
provides nearly continuous space observations of solar oscillation modes with angular degree 
t ranging from to ~ 300. The medium-^ data represent Doppler velocity images acquired 
with 1-minute cadence and a spatial sampling of 10" (0.6 heliographic degrees per pixel) after 
some averaging onboard the spacecraft and ground calibration. These data can be used for 
imaging the solar tachocline and studying its evolution with the solar cycle. 
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We select the MDI observations from July 25 through August 22, 1996, covering the 
whole Carrington rotation CR1912. This period is during the solar minimum between so- 
lar cycles 22 and 23. This helps to avoid potential complications caused by strong surface 
magnetic fields present during other periods of the solar cycle. For each day of the observa- 
tions (1440 minutes long), we remap the data into heliospheric Carrington coordinates using 
Postel's projection. The main purpose of this remapping is to filter out the signals of solar 
convection and /-mode oscillations. For each observation day, we calculate and save the 
cross-correlation functions for both the surface- and deep-focusing measurements. Then, the 
daily cross-correlations are averaged for the whole Carrington rotation, and Gabor wavelet 
fitting is performed. This procedure improves the signal-to-noise ratio. The travel-time 
inversions are done for this Carrington rotation for the surface- and deep-focusing schemes 
separately, and their combination, after the corresponding reference travel times computed 
from the unperturbed simulation (shown in Fig. [2]) are subtracted. 



4.2. Measurement and Inversion Results 

The measured acoustic travel times and inversion results from the surface- focusing 
scheme, averaged over the whole Carrington rotation, are presented in Fig. [9j The 2D 
travel time image has many thin vertical stripes because of noise in the reference travel 
times calculated from the simulated Sun data. A significant travel time dip can be easily 
identified at the annulus radius between ~ 45° and ~ 60°, corresponding to the sound-speed 
perturbation bump near the tachocline region (Fig. [9b). The 2D image of the sound-speed 
perturbations obtained by inversion (Fig. [9b) clearly shows a bump at the location of ap- 
proximately 0.67Rq. We find that the perturbation at the base of the convection zone is not 
latitudinally uniform, being weaker near low latitudes and stronger in latitudes higher than 
about 20°. 

The results for the deep-focusing scheme are presented in Fig. [TUJ The noise level is 
significantly larger for this scheme than for the surface-focusing scheme. Therefore, we ap- 
plied a 3-pixel boxcar average to the 2D travel-time image (Fig. [TUk). and used the averaged 
data for inversion. Similar to the surface-focusing results, the acoustic travel times reveal 
a significant dip for annulus radii from about 46° to 70°. Besides this dark vertical stripe, 
one can also identify in Fig. [TUk white and dark curved band features in both hemispheres, 
extending from an annulus radius of 6° to nearly 40° and between 30° and 60° in latitude. 
It is not clear to us whether these features are artifacts or caused by sound-speed pertur- 
bations in the Sun. The inversion sound-speed results (Fig. [TUb-d) are quite similar to the 
surface-focusing scheme inferences. The sound-speed bump can be clearly identified, and the 
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perturbation exhibits an apparent latitudinal dependence as well. 

The inversion results from a combination of both surface- and deep-focusing measure- 
ments are shown in Fig. [TTJ Because the surface-focusing measurements have a better signal- 
to-noise ratio and cover higher latitudes than the deep-focusing scheme, we intentionally give 
a higher weighting to the surface-focusing measurements in these inversions. The selection 
of this weighting parameter is quite arbitrary, and the weighting has been chosen as 60% 
for the surface- focusing and 40% for the deep-focusing measurement. It seems that such 
combined measurements are useful, but may still need further optimization. 

It is important to compare our sound-speed profiles with the results from global he- 
lioseismology obtained by inversion of normal mode frequencies. Figure [12] shows the com- 
paris qn of the global heliosei smology result obtained from the same SOHO/MDI medium-^ 



data (IKosovichev et al.lll997l ) with the sound-speed perturbation inferred from the combined 
measurements in this study. It is found that both results are in good agreement in terms 
of location and magnitude of the sound-speed bump near the tachocline, as well as in the 
convection zone. 

The good agreement of the latitudinally averaged sound-speed perturbation with the 
previous global helioseismology result indicates that our measurements are generally correct. 
However, the latitudinal variation of the perturbation is surprising . It certainly requires 



further investigation, in particular, of MDI instrumental distortion (IKorzennik et al.ll2004l ) 
and non-uniform MTF effects. 



5. Discussion 

Using numerically simulated global acoustic wavefields and solar oscillation observations 
from SOHO/MDI, we have demonstrated for the first time that time-distance helioseismology 
is capable of measuring sound-speed perturbations in the tachocline located about 200 Mm 
beneath the solar photosphere that are relatively small in magnitude and size. This local 
helioseismology method provides 2D images of these perturbations as a function of latitude 
and depth. This newly developed technique has an advantage over global helioseismology in 
that it can distinguish sound-speed perturbations of different amplitude located at different 
latitude, and is not limited to the North-South symmetrical part. Although we only present 
in this paper an example with maximum perturbations of 0.7% and 0.5%, our experiments 
demonstrate that the technique is capable of measuring the North-South difference in the 
sound-speed structure of the tachocline as small as 0.05%. 

Our inversions of the numerical simulation data have shown that although we can recover 
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the major features in the tachocline region, the inferred sound-speed perturbations are more 
widely spread into areas beneath the tachocline. We think this may be caused by the ray- 
approximation sensitivity kernels used in our inversions or by the relatively high realization 
noise in the simulated wavefield. We plan to continue the work on improving the resolution 
of the deep interior beneath the tachocline. 

Our initial analysis of the MDI data shows a latitudinal dependence of sound-speed 
perturbation at the tachocline. We believe this needs further investigation and modeling 
of potential instrumental effects. Furthermore, it will be extremely useful to investigate 
how the tachocline varies with the solar cycle, and we plan to continue these studies in this 
direction. 
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Fig. 1. — Sketches of wave paths in a meridional plane for two measurement schemes 
surface-focusing and b) deep-focusing. 
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Fig. 2. — Acoustic phase travel time measured from the numerically simulated global wave- 
fields using solar model S. The distance step is 0?6 in the heliographic coordinate. 
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Fig. 3. — a) Sound-speed perturbation model in the tachocline region used in the numerical 
simulation of global acoustic wavefields. The perturbation maxima, 0.7% and 0.5%, are 
centered at O.7O_R and latitude of 30° in Northern and Southern hemispheres, respectively. 
The white dashed curve indicates the radius of 0.70i? Q (the bottom of the convection zone), 
b) Latitudinally averaged sound-speed perturbation of this tachocline model, with the black 
curve from the Northern hemisphere, and the red from the Southern. 



- 15 - 



a) 



Travel Time (min) 



T3 
3 



to 




26 46 66 
Annulus Radius (deg) 



6 26 46 66 86 
Annulus Radius (deg) 



c) 




d) 



15N 



o 

'O 



15S 



0.006 



0.004 



0.002 



0.000 



-0.002 




0.4 0.5 0.6 0.7 0.8 0.9 1.0 
r/R 



Fig. 4. — Measurements and inversion results for the surface-focusing measurement scheme: 
a) acoustic travel times measured from the numerical simulation with the tachocline model, 
relative to the travel times without tachocline; b) latitudinally averaged relative travel times, 
with dark from the Northern hemisphere and grey from the Southern; c) 2D inversion results 
for the sound-speed perturbation; d) latitudinally averaged inversion results, with solid lines 
showing inversions and dotted lines showing the averaged perturbation model. Dark lines 
are for the Northern hemisphere and the red lines are for the Southern. Errors bars indicate 
2a formal error estimates. 
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Fig. 5. — A sample of the averaging kernels for the surface- focusing inversions for the target 
locations at 18° latitude and a) 0.70i? Q , b) 0.80i? Q , and c) 0.90i? Q . The cross in each panel 
indicates the targeted location. 
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Fig. 6. — Same as Fig. HI but for the deep-focusing scheme. Note that this scheme provides 
less coverage than the surface-focusing scheme. 
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Fig. 7. — Same as Fig. but for the deep-focusing scheme. 
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Fig. 8. — Same as panels c) and d) in Fig. HI but for the combined inversion of the surface- 
and deep-focusing measurements. 
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Fig. 9. — Same as Fig. HI but for the SOHO/MDI medium-^ solar oscillation data and 
without showing separate hemispheric results. 
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Fig. 10. — Same as Fig. [9J but for the deep-focusing scheme. 
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Fig. 11. — Same as panels c) and d) in Fig. [9j but for the combined inversion of the surface- 
and deep-focusing measurements. 
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Fig. 12. — Comparison of th e sound-speed perturba tion obtained from global helioseismol- 
ogy (symbol with error bars; iKosovichev et al.lll997l ). and the inversion results obtained by 
combining surface- and deep- focusing measurements of acoustic travel times (solid curve). 



